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In conventional lasers, the optical cavity that confines the photons 
also determines essential characteristics of the lasing modes such as 
wavelength, emission pattern, directivity, and polarization. In random 
lasers, which do not have mirrors or a well-defined cavity, light is 
confined within the gain medium by means of multiple scattering. 
The sharp peaks in the emission spectra of semiconductor powders, 
first observed in 1999, has therefore lead to an intense debate about 
the nature of the lasing modes in these so-called lasers with resonant 
feedback. In this paper, we review numerical and theoretical studies 
aimed at clarifying the nature of the lasing modes in disordered 
scattering systems with gain. The last decade has witnessed the 
emergence of the idea that even the low-Q resonances of such open 
systems could play a role similar to the cavity modes of a conven- 
tional laser and produce sharp lasing peaks. We will focus here 
on the near-threshold single-mode lasing regime where nonlinear 
effects associated with gain saturation and mode competition can 
be neglected. We will discuss in particular the link between random 
laser modes neai^ threshold and the resonances or quasi-bound (QB) 
states of the passive system without gain. For random lasers in 
the localized (strong scattering) regime, QB states and threshold 
lasing modes (TLM) were found to be nearly identical within the 
scattering medium. These studies were later extended to the case of 
more lossy systems such as random systems in the diffusive regime 
where it was observed that increasing the openness of such systems 
eventually resulted in measurable and increasing differences between 
quasi-bound states and lasing modes. Very recently, a theory able 
to treat lasers with arbitrarily complex and open cavities such as 
random lasers established that the threshold lasing modes are in fact 
distinct from QB states of the passive system and are better described 
in terms of a new class of states, the so-called constant-flux states. 
The correspondence between QB states and lasing modes is found 
to improve in the strong scattering limit, confirming the validity of 
initial work in the strong scattering limit. © 2010 Optical Society of 
America 
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1. Introduction 



The investigation of laser action in complex media with disorder has a long history go- 
ing back to the early days of laser physics (for a review, see [1-5]). Beginning in the mid 
1990's there was a resurgence of interest in this topic both for its intrinsic interest and 
because of a possible relation to the phenomenon of Anderson localization [6] previ- 
ously studied mainly in the context of electronic systems. Random lasers are disordered 
media with gain that do not possess a light-trapping cavity beyond the confinement pro- 
vided by multiple-scattering from the disorder itself. Hence they are usually extremely 
open, low finesse lasers. Initially it was unclear whether such systems could produce 
naiTow lasing lines without any well-confined electromagnetic modes, and while initial 
experimental studies did find strong amplification near the transition frequency deter- 
mined by the gain medium, discrete lines were not observed [7-9]. Subsequent studies 
in smaller systems with focused pumping did find discrete lasing lines not necessar- 
ily located at the center of the gain curve and photon statistics characteristic of gain 
saturation [10-14], demonstrating that in some cases random lasers (RLs) behave very 
much like conventional multimode lasers except for their relatively high thresholds and 
their pseudo-random emission patterns in space. The experimental observations of laser 
peaks have naturally called for the search for a feedback mechanism leading to light 
trapping within the scattering medium. There is in fact a case where light can be well- 
confined inside an open disordered medium. Such confinement occurs when the scatte- 
ring is extremely strong and the system is in the regime of Anderson localization [15]. 
However, except in quasi- ID geometries [16], the vast majority of experiments on RLs 
do not appear to be in the localized regime, so the question of whether laser action in a 
diffusive (L ^> tj or quasi-ballistic (L ~ tj medium has a qualitatively different nature 
and origin with respect to conventional lasers remained open for some time (here L is 



the system size and I is the optical elastic mean free path). 

With the renewed experimental interest in RLs came also a number of attempts to 
generalize laser theory to describe such a system. Early on a major distinction was made 
between conventional lasers which operate on resonant feedback and RLs which at least 
in some cases were supposed to operate only on non-resonant feedback (NRF) [4]. In 
the case of NRF the light intensity in the laser was described by a diffusion equation 
with gain but the phase of the light field and hence interference did not play a role. A 
key finding is that there is a threshold for amplification when the diffusion length for 
escape ~ dlJ jl became longer than the gain length (here J = 2,3 is the dimen- 
sionality). The spatial distribution of intensity above threshold would be given by the 
solution of a diffusion equation. In this approach there would be no frequency selec- 
tivity and the amplified light would be peaked at the gain center. Clearly such a de- 
scription would be inadequate to describe RLs based on Anderson localized modes, as 
such modes ai^e localized in space precisely due to destructive interference of diffusing 
waves arising from multiple scattering. 

This question itself is related to a basic question in non-linear optics: can a system, 
disordered or not, which is so leaky that it has no isolated linear scattering resonances 
nonetheless have sharp laser lines due to the addition of gain? And if so, how are the 
modes associated with these lines related to the broad and difficult to observe reso- 
nances of the passive cavity? For an open diffusive or quasi-ballistic medium in two 
or three dimensions the resonance spacing in wavevector will decrease as X'^^^ /L'^, 
whereas the linewidth will scale as l/L^ (diffusive) or as 1/L (ballistic). Therefore 
(unless £k, X ind = 2) the disordered passive cavity resonances strongly overlap and 
cannot be directly observed in linear scattering. 

In the search for a feedback mechanism responsible for the sharp laser peaks ob- 
served experimentally [17], different scenarios have been proposed. As an alternative 
to the early picture of closed scattering loops, the probability of having ring-shaped 
resonators with index of refraction larger than average in the diffusive regime was 
calculated and shown to be substantially increased by disorder correlation due to finite- 
size scatterers [18, 19]. Another scenario was put forward where spontaneously emitted 
photons accumulate gain along very long trajectories. This follows the observation of 
random spikes in the emission spectrum of weakly active scattering systems in single 
shot experiments [20,21]. These "lucky photons" accumulate enough gain to activate a 
new lasing mode with a different wavelength after each excitation shot. The experimen- 
tal study of the modal decay rates in microwave experiments leading to the observation 
of anomalous diffusion has brought forward the existence of longer-lived prelocalized 
modes in an otherwise diffusive system [22]. An experimental indication of the coexis- 
tence of extended and locaUzed lasing modes has been presented recently [23]. It was 
suggested that these longer lived-modes could be responsible for lasing. However, al- 
though they are possibly achieved in some specific situations, those different scenarios 
cannot explain the whole set of experimental observations 

In this paper, we present recent work, both numerical and analytical, which has 
shown that within semiclassical laser theory, in which the effects of quantum noise 
are neglected, definite answers to these questions can be given, without resorting to ex- 



otic scenaiios. Sharp laser lines based on interference (coherent feedback) do exist not 
only in strongly scattering random lasers where the localized regime is reached [24-26] 
but also in diffusive random lasers (DRLs) [27,28], and even for weak scattering [29]. 
Numerical studies have shown that they ai^e associated with thi^eshold lasing modes 
(TLMs), which, inside the cavity, are similai" to the resonances or quasi-bound (QB) 
states of the passive system (also called quasi-normal modes). The resemblance is ex- 
cellent in the localized case [25,26] and deteriorates as scattering is reduced. A new 
theoretical approach based on a reformulation of the Maxwell-Bloch equations to ac- 
cess the steady state properties of arbitrarily complex and open cavities allows one to 
calculate the lasing modes in diffusive and even in weakly scattering random lasers 
{£ ~ L) [30-34]. A major outcome of this approach is the demonstration that although 
lasing modes and passive modes can be very alike in random systems with moder- 
ate openness in agreement with the above numerical results, they feature fundamental 
differences. Their distinctness increases with the openness of the random system and 
becomes substantial for weakly scattering systems. Constant-flux (CF) states are in- 
troduced which better describe TLMs both inside and outside the scattering medium 
for any scattering strength. In addition this theoretical approach allows one to study 
the multimode regime in DRLs and get detailed information about the effects of mode 
competition through spatial hole-burning, which appear to differ from conventional 
lasers. 

In this last decade, different types of random lasers (semiconductor powders, pi- 
conjugated polymers, scattering suspension in dyes, random microcavities, dye-doped 
nematic liquid crystals, random fiber lasers, ...) have been considered in the literature. 
We will focus throughout this review mostly on two dimensional (2D) random lasers 
which consist of randomly distributed dielectric nanoparticles as scatterers. This choice 
makes possible the numerical and theoretical exploration of 2D finite-sized opened 
samples where transport can be ballistic, diffusive (in contrast to ID) or localized [35], 
by adjusting the index contrast between the scatterers and the background medium. 

The outline of this review is as follows: In section 2 we review early numerical 
explorations of localized and diffusive random lasing demonstrating the existence of 
threshold lasing modes (TLMs) in all regimes. In section 3 we present recent numeri- 
cal work based on a time-independent model, which indicates the difference between 
passive cavity resonances and threshold lasing modes, discussing only single-mode 
random lasing. The following section will explain why in principle, quasi-bound (QB) 
states cannot describe TLMs. The last section will introduce the concept of constant- 
flux (CF) states and describe the self-consistent time-independent approach to describe 
random lasing modes at threshold as well as in the multimode regime. 

2. Early numerical explorations: Time-dependent model 

2.1. Localized case 

From a modal point of view, Anderson localization means that for strong disorder, the 
eigenmodes of the wave equation are spatially localized in a volume of finite size 2^, 
where ^ is the localization length. More precisely, they are spatially localized solutions 



of the Maxwell equations with tails, which decay exponentially from their center, 
being the decay length. In the case of scattering particles, the value of the localization 
length is controlled by the index contrast between the particles and the background 
medium, the size of the particles, the optical wavelength and the amount of disorder. 
In practice, when finite-size systems in the localization regime are considered, two 
opposite cases may occur : (1) > L and (2) < L where L is the system size. In 
the first case, the system is not lai^ge enough for the light to be confined by disorder 
within the volume of the system. In case (2), light is localized since it cannot escape 
domains larger than t, . More precisely, localized modes are coupled to the boundaries 
via their exponential tails. The leakage rate of an exponentially localized quasi-bound 
(QB) state varies as exp(— 2r/(^) with r the distance to the boundaries [36] . Hence, in 
sufficiently large systems QB states located far from the boundaries (which constitutes 
majority of the QB states except a fraction proportional to t, /L) feature a very small 
leakage i.e. a good quality factor. 

In this subsection, we will consider case (2). Localized modes in a disordered scat- 
tering system are quite alike the modes of standard optical cavities, such as the Fabry- 
Perot. [37]. Hence, one can expect that in the presence of gain, the lasing modes in 
this regime of strong disorder will be close to the localized QB states of the passive 
system without gain, in the same way as the lasing modes of a conventional cavity are 
built with the QB states of the passive cavity. In order to verify that this is really the 
case one must have access to the individual modes of both the passive system and the 
active system. Experimentally, such a demonstration has not been achieved yet, essen- 
tially because the regime of Anderson localization is difficult to reach and to observe 
in optics. Besides, until recently there was no fully developed theory describing ran- 
dom lasing modes and their relationship with the eigenstates of the passive system. The 
easiest way to check this conjecture has been to resort to numerical simulations. 

Historically, most of the early numerical studies of random lasers were based on 
the diffusion equation (see references in [4]). However, it is not possible to take into 
account under the diffusion approximation the interference phenomena which are at 
the heart of Anderson localization. This is why Jiang and Soukoulis [24] proposed to 
solve the time-dependent Maxwell equations coupled with the population equations of 
a four-level system [38]. The populations Ni,i = 1 to 4 satisfy the following equations 



where Wp is the rate of an external mechanism which pumps electrons from the fun- 
damental level (1) to the upper level (4). The electrons in level 4 relax quickly with 
time constant T43 to level 3. The laser transition occurs from level 3 to level 2 at fre- 
quency (Oa- Hence, electrons in level 3 can jump to level 2 either spontaneously with 
time constant T32 or through stimulated emission with the rate (E/hcOa)dP/dt. E and P 
are the electric field and the polarization density respectively. Eventually, electrons in 



dNi/dt =N2/i:2i-WpNi 

dN2/dt =N3/T32-N2/T21 - / h(Oa)dV / dt 

dN^/dt = Na/za3 -N3/T32 + {E/n(Oa)dV/dt 
dN^/dt = -A^3/t43 + WpiVi 



(1) 
(2) 
(3) 
(4) 



level 2 relax quickly with time constant T21 from level 2 to level 1. In these equations, 
the populations A'^,, the electric field E and the polarization density P are functions of 
the position r and the time t. 

The polarization obeys the equation 



where AA/^ = N2—N3 is the population density difference. Amplification takes place 
when the rate Wp of the external pumping mechanism produces inverted population 
difference AA'^ < 0. The line width of the atomic transition is AcOa = I/T32 + 2/72 where 
the collision time T2 is usually much smaller than the lifetime T32. The constant K is 
given by K- = 3c^/2(0^t:32 [38]. 

Finally, the polarization is a source term in the Maxwell equations. 



The randomness of the system arises from the dielectric constant £(r), which depends 
on the position r. This time-dependent model has been used in random ID systems 
consisting of a random stack of dielectric layers separated by gain media [24] and in 
random 2D systems consisting of a random collection of circular particles embedded 
in a gain medium (Fig. 1) [25]. In both cases, a large optical index contrast has been as- 
signed between the scatterers and the background medium to make sure that the regime 
of Anderson localization was reached. The Maxwell equations are solved using the 
finite-difference time-domain method (FDTD) [39]. To simulate an open system, per- 
fectly matched layers are introduced at the boundaries of the system [40]. The pumping 
rate Wp is adjusted just above lasing threshold in order to remain in the single-mode 
regime. 

In ID, the QB states of the passive system were obtained independently using a time- 
independent transfer matrix method [41]. In 2D, the Maxwell equations were solved 
without the polarization term in (7) using again the FDTD method. First, the spectrum 
of eigenfrequencies was obtained by Fourier transform of the impulse response of the 
system. Next, QB states were excited individually by a monochromatic source at each 
of the eigenfrequencies. 

Finally, in ID systems [41] as well as in 2D systems [26], lasing modes obtained by 
the full time-dependent model with gain and localized QB states of the corresponding 
passive system without gain were compared and found to be identical with a good pre- 
cision. This was verified for all modes obtained by changing the disorder configuration. 
An example of a 2D lasing mode and the corresponding QB state of the same system 
(Fig. 1) without gain are displayed in Fig. 2. These results confirmed that the QB states 
of a localized system play a role similar to the eigenmodes of the cavity of a conven- 
tional laser. The only difference is the complicated and system-dependent nature of the 
localized modes as opposed to the well-known modes of a conventional cavity. These 
results are in good agreement with the theoretical results described in Section 5, which 
show that inside systems in the localized regime, the single lasing modes just above 
threshold are close to the high-Q resonances of the passive system. 



d^P/dt^ + A(OadP/dt + (of,P = k.AN.E 



(5) 



dH/dt = -cV X E 
e{r)dEdt = cV X H - 4nd'P/dt. 



(6) 
(7) 



FIGURE 1. (Color online). Example of a random realization of 896 circular scat- 
terers contained in a square box of size L=5 microns and optical index n=l. The 
radius and the optical index of the scatterers are respectively r=60 nm and n=2 . 
The total system of size 9 microns is bounded by perfectly matched layers (not 
shown) in order to simulate an open system 

2.2. Diffusive case 

We have seen in the previous section that random lasers in the Anderson localization 
regime should behave like conventional lasers. They should exhibit discrete laser peaks 
above threshold in agreement with the experimental observations of laser action with 
resonant feedback. However, subsequent measurements of the mean free path showed 
that none of the experimental cases that displayed discrete laser peaks were in the 
localized regime. Instead, they were found to be in the diffusive regime and some even 
in the quasi-ballistic regime [29]. In such systems, there are no localized modes so that 
the observation of laser action with resonant feedback has been the subject of much 
debate. 

Only very recently, numerical evidence was given that even diffusive systems with 
low Q resonances could exhibit lasing with resonant feedback [27]. The random 2D 
systems described in the previous section consisting of random collections of circu- 
lar particles embedded in a gain medium have been investigated with the same time- 
dependent model. To be in the diffusive regime instead of the localized regime, a 
smaller optical index contrast 5„ = 0.25 instead of 5„ = 1.0 has been assigned between 
the scatterers and the background medium. Solving the Maxwell equations coupled to 
the population equations, laser action characterized by a sharp peak in the emission 
spectrum was observed just above a threshold, albeit high. An example of the corre- 
sponding lasing mode is displayed in Fig. 3a. In contrast to the localized case, the las- 



FIGURE 2. (Color online), (a) Spatial distribution of the amplitude of a lasing 
mode in the localized regime (« = 2) and (b) that of the corresponding QB states of 
the same random system without gain. The square delimits the scattering medium. 
The amplitude rather than the intensity is represented for a better display of the 
small values of the field. 

ing mode is now extended over the whole system. Moreover this is a complex mode in 
the sense that it contains a substantial traveling wave component [27]. However, in this 
work comparison of the lasing modes with the QB states of the passive cavity could not 
be caiTied out by using the time domain method as it was done in the localized regime. 
Due to strong leakage through the boundaries, resonances are strongly overlapping in 
the frequency domain and one cannot excite them individually by a monochromatic 
source. 

To circumvent this difficulty, an indirect method has been used to compare the lasing 
modes with the resonances of the passive system. This method is inspired by the Fox-Li 
modes, which in conventional laser physics are modes of an open cavity [38,42,43]. 
The Fox-Li modes are field distributions whose profile is self-repeating in a complete 
round trip within the Fabry-Perot laser cavity while decaying because of the diffraction 
losses due to finite surface area of the end miiTors. Analogously, if the lasing modes of 
the diffuse system are related to the resonances of the passive system, they should decay 
by self-repeating themselves when pumping and population inversion are turned off. 
To study the evolution of the mode profile with time, the following spatial correlation 
function was introduced [27] 



which compares the mode profile S'{r,t) at time t with the mode profile at the 
initial time tQ. Here, ^ is the scattering medium. The field has been normalized 
(^(j,t) = £'(r,f)/[/ /^J^r£'^(r,f)]^/^ to counterbalance the decay due to the leakage 
through the boundaries. This correlation function oscillates at the laser frequency be- 




(8) 




FIGURE 3. (Color online), (a) Spatial distribution of the amplitude of a lasing 
mode in the diffusive regime, (b) Spatial distribution of the field amplitude after 
the pump has been stopped and the polarization term has been set to zero. The 
spatial distribution of scatterers is the collection shown in Fig. 1 but here, the 
optical index of the scatterers is n = 1 .25 instead of n = 2 in Fig. 2 

tween -1 and +1 if the normalized mode profile is recovered at each period (Fig. 4). 
Otherwise, the amplitude of the oscillations should decay with time. This correlation 
function was used in [27] to check whether the first lasing mode at threshold for dif- 
fusive random laser indeed corresponds to a Fox-Li mode of the passive system. The 
pumping is set to zero after the lasing mode has been established so that at later times 
the field can evolve by itself. The long time evolution of the spatial correlation function 
associated with this free field is displayed in Fig. 5a. The decay of the total energy of 
the system is also shown. While energy decay is observed over 6 orders of magnitude, 
the spatial correlation function is seen to oscillate between values close to -1 and -i-l 
meaning that the initial lasing mode profile Sir, to) is reproduced at each period with a 
good accuracy. The decaying field amplitude has the spatial distribution which is shown 
in Fig. 3b until eventually, the correlation function decays to zero when the decaying 
field reaches the noise level. This result demonstrates that the TLM is very close to 
a resonance of the passive system, when measured inside the scattering medium. For 
comparison, the evolution of the spatial coiTclation function for an initial field created 
by an arbitrary distribution of monochromatic sources at the laser frequency is dis- 
played in Fig. 5b. The fast decay of Cg{tQ,t) after the sources have been turned off 
indicates that this field distribution is not a QB state of the passive system. 

The decay rate observed coiTcsponds to a quality factor of 30, to be compared with 
the value 10^ found in the localized case. This result shows that a "bad" resonance in 
a leaky disordered system can nevertheless turn into a lasing mode in the presence of 
an active medium. This result is in stark contrast with the common belief that random 
lasing with resonant feedback involves the presence of resonances with high quality 
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FIGURE 4. (Color online). Short time behavior over a few cycles of the correla- 
tion function, C^(fo,f)' for W ^ localized lasing mode as in Fig. 2 and for (b) a 
diffusive lasing mode as in Fig. 3. The periodic square function in (a) is typical of 
a standing wave while the sinus-like function in (b) is characteristics of a traveling 
wave. [27] 



factors. It provides a consistent explanation for the experimental observation of random 
lasing with resonant feedback even far from the localized regime, without resorting to 
other scenarios such as those reviewed in the introduction [17-21]. 

The comparison of patterns between Fig. 3a and Fig. 3b shows that the lasing mode 
and the QB modes are close to each other inside the scattering system as confirmed 
by the evolution of the correlation function, which has been defined only inside the 
system. However, one also notices that outside the scattering medium, the field dis- 
tributions differ substantially. The free propagating field outside the scattering system 
in Fig. 3b reproduces the laser field distribution in Fig. 3a with significant distortions 
due to the enhancement of the amplitude towards the external boundaries of the total 
system. Hence, the comparison between both figures indicates that if the lasing modes 
and the QB modes are similar inside the scattering system, they differ noticeably out- 
side. Isloreover, a careful examination of the correlation function in Fig. 5a shows that 
it oscillates between two extremal values, which slowly depart from -1 and -i-l well 
prior to the ultimate fast decay. This is in contrast with the long time behavior of the 
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FIGURE 5. (Color online). Correlation function (full line) and energy decay 
(dashed line) vs time of (a) the lasing mode when the pump is turned off and 
(b) an arbitrary field distribution at the frequency of the lasing mode 



correlation function in the localized regime (not shown), which displays oscillations 
between -1 and +1 with a very good precision for time scales much longer than the 
time scale in Fig. 5a. This result indicates also that inside the scattering system, the 
lasing mode is close to but not identical to a QB state. 

In conclusion, the time dependent model has provided direct evidence of the close- 
ness of lasing modes and passive cavity resonances at least in the localized case. In the 
diffusive regime, the lasing modes are also found rather close to the QB modes although 
small discrepancies manifest themselves. We also found that this holds inside the scat- 
tering medium. Outside the scattering system however, differences become more sig- 
nificant. The advantage of the time dependent model is that one has access in principle 
to the full non-linear dynamics of the laser system. However QB states with low qual- 
ity factors are not accessible with this approach. Hence, the measure of the difference 
between TLM and QB states has been indirectly achieved by using the spatial correla- 
tion function. Another limitation of this method is related to the various time constants 
involved in this model which lead to time consuming computations, particularly when 
one wishes to vary disorder and study an ensemble of disorder configurations. To over- 
come these limitations, different approaches such as solving the wave equation in the 



frequency domain have been used. Several approaches of this kind will be described in 
the next section [13,44-46]. The recent theoretical approach based on a different class 
of states, the so-called constant flux (CF) states, and taking into account non-linear 
interactions will be described in section 5. 

3. Numerical simulations: Time-independent models 

Different models have been proposed in the frequency domain to solve the wave equa- 
tion. In ID, it is possible to employ the transfer matrix method similar to that used 
in [41] for studying the lasing modes in an active layered random system. A direct 
comparison between TLMs and quasi-bound states of the corresponding passive ran- 
dom system is proposed in the first part of this section. In 2D, the multipole method 
has been used, which also provides a direct comparison of the QB states and the las- 
ing modes of a 2D-disordered open system. The comparison presented in the second 
part of this section has been carried out for refractive index of the scatterers ranging 
from n'l = 2.0, (localized regime) to ?ij = 1.25 (diffusive regime). We alternatively used 
a different approach based on the finite element method to obtain the passive modes, 
which turned out to be much less computationally demanding in the weakly scattering 
regime. A brief description of both methods is provided in Appendices A and B. 

3. 1 . One-dimensional random lasers 

Employing the transfer matrix method, similar to that used in [41], we study the lasing 
modes in a one-dimensional (ID) random system and compare them with the QB states 
of the passive random system. The random system is composed of 161 layers. Dielec- 
tric material with index of refraction «i = 1.05 separated by air gaps («2 = 1) results in 
a spatially modulated index of refraction n{x). Outside the random medium no = I. The 
system is randomized by specifying thicknesses for each layer as 2 = (<ii,2) ( 1 + 0' 
where {d\) = 100 nm and {di) = 200 nm are the average thicknesses of the layers, 
T] = 0.9 represents the degree of randomness, and is a random number in (-1,1). The 
length of the random structure L is normalized to (L) =24100 nm. Linear gain is sim- 
ulated by appending an imaginary part to the dielectric function e{x) = e'{x) + ie"{x), 
where s'{x) = n^{x). This approximation is valid at or below threshold [46]. The com- 
plex index of refraction is given by h{x) = y£{x) = n'{x) + in" , where n" < 0. We 
consider n" to be constant everywhere within the random system. This yields a gain 
length Ig = \ l/k"\ = \/\n"\k {k = 2n/X is the vacuum frequency of a lasing mode) 
which is the same in the dielectric layers and the air gaps. The real part of the index of 



We find the frequency k and threshold gain k of each lasing mode within the wave- 
length range 500 nm < A < 750 nm. The results are shown in Fig. 6. Finding 'match- 
ing' QB states for lasing modes with large thi^esholds (large \k |) is challenging due to 
large shifts of the solution locations [Fig. 6(region c)]. However, there is a clear one- 
to-one correspondence with QB states for the lasing modes remaining [Figs. 6(region 
a) and (region b)]. It is straightforward to find the matching QB states for these las- 
ing modes and calculate their differences. The average percent difference between QB 



refraction is modified by the imaginary part as n'{x) 
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FIGURE 6. The frequencies k of quasi-bound modes (crosses) and lasing modes 
with linear gain (open diamonds) together with the decay rates of QB states 
and the lasing thresholds k" of lasing modes. The horizontal dashed lines separate 
three different regions of behavior: (a) lasing modes are easily matched to QB 
states, (b) clear differences appear but matching lasing modes to QB states is still 
possible, (c) lasing modes have shifted so much it is difficult to match them to 
QB states. The QB state with the largest decay rate and the lasing mode with the 
largest threshold are circled, though they may not be a matching pair 



State frequencies and lasing mode frequencies in Fig. 6(region a) is 0.013% while it 
is 0.15% in Fig. 6(region b). The average percent difference between QB state decay 
rates and lasing thresholds k in Fig. 6(region a) is 2.5% and in Fig. 6(region b) is 
21%. 

The normalized intensities of the QB states Iqb and lasing modes with linear gain 
Ilg are also compared. Figure 7 shows representative 'pairs' of modes from the three 
regions shown in Fig. 6. The spatially averaged relative difference between each pair 
of modes is calculated by 

For small thresholds [Fig. 7(a)], the difference between the lasing modes and the 
matching QB states is very small. The average percent difference between all pairs of 
modes in this region is (a^) = 4.2%. For lasing modes with slightly larger thresholds 
[Fig. 7(b)], there are clear differences. Nevertheless, we may confidently match each 
lasing mode in this region with its corresponding QB state. The average percent differ- 
ence between all pairs of modes in this region is (aj) = 24%. As mentioned earlier, it 
is challenging to find matching pairs of lasing modes and QB states for large thresh- 
olds. Figure 7(c) compares the lasing mode with the largest threshold and the QB state 
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FIGURE 7. (Color online) Spatial intensity distributions of quasi-bound modes 
Iqb{x) (red solid lines) and lasing modes Ilg{x) (black dashed lines) from each of 
the three regions in Fig. 6. Representative samples were chosen for each case, (a) 
The lasing mode intensity is nearly identical to the QB state intensity with Od = 
1.7%. (b) A clear difference appears between the lasing mode and the QB state, 
with = 21.8%, but they are still similar (c) The lasing mode with the largest 
threshold and QB state with the largest decay rate are compared, with Od = 198%. 
Though these two modes are fairly close to each other [circled in Fig. 6(region 
c)], their intensity distributions are quite different. 

with the largest decay rate [circled in Fig. 6(region c)]. Though these two modes are 
fairly close to each other in terms of k, k^, and k , their intensity distributions are quite 
different. Indeed, there may be no coiTcspondence between the two. 

The deviation of the lasing modes from the QB states can be explained by the mod- 
ification of the transfer matrix. In the passive system, k^ is constant, but k i = kQn{x) 
varies spatially. With the introduction of gain, k becomes constant within the random 
system, and feedback due to the inhomogeneity of k is removed. However, introducing 
gain generates additional feedback inside the random system caused by the modifica- 
tion in the real part of the wave vector k' = kn'{x). Neglecting this effect results in some 
coiTcspondence between lasing modes and QB states even at large thresholds [47]. Fur- 
thermore, since there is no gain outside the random system, k suddenly drops to zero 
at the system boundary. This discontinuity of k generates additional feedback for the 
lasing modes. In this weakly scattering system, the threshold gain is high. The large 
drop of k at the system boundary makes the additional feedback stronger. 

3.2. Two-dimensional random lasers 

We turn now to the 2D case. A different approach based on the multipole method 
has been used. The multipole method is best suited to characterize multiple scattering 



problems involving scatterers with circulai^ cross-section. This method has been used 
to compute the scattering of a plane wave by a random collection of cylinders [29,48], 
to calculate the defect states in photonic crystals [49], to construct the exact Green's 
function of a finite system [50], or to calculate the local density of states [51]. This 
method has also been used to explain the anomalously large Lamb shift that occurs 
in photonic crystals by calculating the QB states in such structure [52]. Finally, The 
multipole method can be used to characterize the modes of three dimensional structures 
composed of cylinders [53] and in particular to find the modes of the phothonic crystal 
fibers [54-56]. It will be used here to calculate the QB states and the lasing modes of 
the 2D-disordered scattering systems of the kind shown in Fig. 1 and studied in the 
previous section for different regimes of scattering. Details about this method can be 
found in Appendix A. 

This method is based essentially on a search for the poles of a scattering matrix. 
Because the system is open, the problem is not Hermitian and hence there are no modes 
occurring for real wavelengths. The poles of the QB states all occur in the complex 
plane at wavelengths X = X' + iX", with causality requiring that X" > 0. The real part 
of the wavelength X' determines the resonance wavelength of the QB state, while the 
imaginary part X" determine the quality factor Q = X' / {2X") of the mode [52]. 

The same method is used to find the lasing modes (TLM) at threshold. It is necessary 
to find this time the poles of the scattering matrix in the two-dimensional space (A', e^') 
of real wavelengths (A" = 0) and imaginary component of the complex dielectric con- 
stant outside the scatterers where the gain is distributed. It can also be used to find the 
lasing modes when gain is localized inside of the scatterers. In this case the poles of 
the scattering matrix are searched in the space of real wavelengths (X" = 0) and the 
imaginary part of the dielectric constant of cylinders e'/) 

The multipole method is both accurate and efficient: the boundary conditions are 
analytically satisfied, thus providing enhanced convergence, particulai^ly when the re- 
fractive index contrast is high. However, in the case of lai^ge systems the method can 
be slow (given that field expansions are global, rather than local) when it is neces- 
sary to locate all poles within a sizable wavelength range. An other extremely efficient 
time-independent numerical method based on Finite Element Method [57] has been 
tested. This method is briefly described in the B. We checked that the results obtained 
by both methods, the (purely numerical) finite element method and the (semi-analytic) 
multipole method were identical with a good precision. 

3.2.a. Localized case 

We first consider the localized case (nj = 2.0) for which a complete comparison of the 
QB states and the lasing modes was possible with the time-dependent FDTD-based 
method (section 2.1), thus providing a reference comparison for the multipole calcu- 
lations. The lasing mode is found at a wavelength X' = 446.335 nm for a value of the 
imaginary part of the refractive index n'j' = —1.967 x 10^"^, representing the pumping 
threshold for this mode. The spatial distribution of its amplitude is shown in Fig. 8(b). 
The QB states of the passive system are calculated in the spectral vicinity of the lasing 
mode. The number of required multipoles was N,nax = 4 (see A). Figure 8(a) shows the 



HGURE 8. (Color online) The intensity \E\^ of the locaHzed QB state (a) 
(Media 1) and corresponding lasing mode (b) (Media 2) calculated using multi- 
pole method for a 2D-disordered scattering system of the kind shown in Fig. 1 
with the refractive index of the cylinders = 2.0. 

QB state which best resembles the lasing mode. Its wavelength and quaUty factor are 
respectively X' = 446.339 nm and Q = 8047. The relative difference between the two 
modes is (aj) = 0.05%. These calculations provide confirmation that the lasing modes 
and the QB states are the same inside the scattering region for high Q-valued states. 




HGURE 9. (Color online). The intensity \E\'^ of the diffusive QB state (a) 
(Media 3) and the lasing mode (b) (Media 4) calculated using the multipole 
method for the same random configuration as in Fig. 8 but with the refractive 
index of the cylinders of «J = 1 .25. 

3.2.b. Diffusive case 

We next consider the diffusive case and choose n\ = 1.25. This is where the time- 
independent method becomes interesting since, in contrast to the FDTD approach, it 
gives a direct access to the QB states. They are accurately calculated in this regime for 
^max = 2 multipoles. Figure 9 shows a lasing mode and its corresponding QB state. 
The lasing mode is found at X' = 455.827 nm for an imaginary part of the refractive 
index = —3.11S x 10^^. The wavelength and the quality factor of the QB state are 
respectively A' = 456.79 nm and Q = 29.2. The lasing mode is therefore red-shifted 
relative to the QB state's wavelength, as a result of the mode-pulling effect. The QB 
state and the lasing mode appear similar in Fig. 9. However, the relative difference 
between the two modes is larger than in the localized case, (a^/) = 14.5%. Fig. 10 
shows the cross-section of the spatial intensity of both modes along x = 2.75. In spite 




of the resemblance, the two profile display visible dissimilarities. This suggests, in the 
diffusive case, that QB states and lasing modes are not exactly the same though they 
exhibit quite similar features. These results are consistent with the findings presented 
in section 2.2. 

3.2.C. Transition case 




FIGURE 1 1 . (Color online) The intensity p of a QB state (a) and a lasing mode 
(b) calculated using multipole method for the same random configuration as above 
but with the refractive index of the cylinders nj = 1 .75. 

It is both informative and interesting to follow the evolution of the lasing modes 
and QB states spatial profile when the index of refraction is decreased progressively, 
allowing to compare systematically the QB states and the random lasing modes (TLM) 
in regime ranging from localized to diffusive. The QB state and lasing modes calcu- 
lated for intermediate cylinder refractive indices n'j = 1.75 and n\ = 1.5 we. displayed 
in Figs. 1 1 and 12. We note that the highly spatially localized mode for n'j = 2 (Fig. 8) 
is replaced for ?ij = 1 .75 by a mode formed by two spatially localized peaks and several 
smaller peaks. For a refractive index of ?ij = 1.5, the mode is still spatially localized 
although on a larger area, but is now formed with a large number of overlapping peaks. 



FIGURE 12. (Color online) The same as in Fig. 1 1 but for = 1 .5. 



A more systematic exploration of the nature of the lasing modes at the transition be- 
tween localized states and extended resonances can be found in [28]. There, a scenario 
for the transition has been proposed based on the existence of necklace states which 
form chains of localized peaks, resulting from the coupling between localized modes. 
The modes shown here support this scenario. It is important to note that the decreasing 
scattering and increasing leakage not only affect the degree of spatial extension of the 
mode but also the nature of the QB states. Indeed, it was shown in [28] that because 
of leakage, extended QB states have a non- vanishing imaginary part associated with a 
progressive component, in contrast to the purely stationary localized states. We present 
in [58] animations of the time oscillation of the real part of the field 9?[*Pexp(— /o)?)] 
of the QB state and of the corresponding TLM for n' = 2 and n' = 1.25. The QB state 
is exponentially decaying in contrast to the lasing mode. The diffusive lasing mode 
clearly exhibit a progressive component, which does not exist in the localized lasing 
mode. 
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TABLE 1. Wavelength X' and quality factor Q of the QB states; lasing frequency 
A' and imaginary part of the refractive index «J' obtained for the threshold lasing 
modes; relative index difference (a^) between QB states and TLM, for four index 
values «' of the scatterers. 

The values of wavelengths and quality factors of the QB states, lasing frequencies 
of the con^esponding threshold lasing modes and associated imaginary part of the re- 
fractive index are summarized in Table 1, together with the relative difference (a^/) as 
defined in Eq.9. 

In order to visualize the increasing difference between TLM and QB states, the cross- 
section of their spatial intensity profile at ;c = 2.75 is plotted in Fig. 13. In Fig. 13(a) 
one cannot distinguish between the lasing mode and the QB state for n' = 1.75, while 
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FIGURE 13. The intensity \E\^ of QB state (blue dashed line) and lasing mode 
(red solid line) at x = 2.75 for «; = 1 .75 (a), n' = 1 .5 (b) and n' = 1 .25 (c). 

for n' = 1.5 (Fig. 13(b)) differences begin to emerge, becoming more pronounced for 
the case of n' = 1.25 (Fig. 13(c)). This is seen also in the increase of the relative dif- 
ference from 5% to 14.5%. Clearly, there is a systematic increase of the discrepancy 
between QB states and lasing modes when index contrast and scattering decrease and 
leakage increases. For very low scattering n' = 1.05, we could not find the QB state 
corresponding to the TLM. Although we may have missed a pole in the complex plane, 
this raises however a serious question on the validity of the comparison of the threshold 
laser mode with QB states when weakly scattering systems are considered. In the next 
section, we will argue that, in principle, QB states cannot be the support of the TLM. 
Section 5 will introduce a different class of states, which offer a valid basis on which 
the TLMs can be described. 

4. Threshold lasing states vs. passive cavity resonances 

Semiclassical laser theory treats classical electromagnetic fields coupled to quantized 
matter and yields the thresholds, frequencies and electric fields of the lasing modes, but 





not their linewidths or noise properties. In order to treat the spatial dependence of lasing 
modes one must go beyond rate equation descriptions and use the coupled non-lineai^ 
Maxwell-Bloch (MB) equations for light coupled to homogeneously broadened two- 
level "atoms" or multilevel generalizations thereof. These equations will be presented 
in section 5 below. While the MB description has been used since the inception of laser 
theory [59,60], in almost all cases simplifications to these equations were made, most 
notably a neglect of the openness of the laser cavity. As random lasers are strongly 
open systems, it is necessary to treat this aspect of the problem correctly to obtain a 
good description of them. 

Historically a first breakthrough in describing Fabry-Perot type lasers with open 
sides was the Fox-Li method [42], which is an integral equation method of finding 
the passive cavity resonances of such a structure. It is widely assumed and stated that 
these resonances or quasi-bound (QB) states are the conect electromagnetic modes of 
a laser, at least at threshold. Often the non-linear laser equations are studied with her- 
mitian cavity modes with phenomenological damping constants representing the cavity 
outcoupling loss obtained, e.g. from a Fox-Li calculation. It is worth noting that there 
are two kinds of cavity loss that occur in lasers; there is the outcoupling loss just men- 
tioned and also the internal absorption of the cavity which can be taken into account 
via the imaginary part of the passive cavity index of refraction. These are very different 
processes as the former describes the usable coherent light energy emitted from the 
laser and the latter simply energy lost, usually as heat, in the laser cavity. 

The QB states of an arbitrary passive cavity described by a linear dielectric function 
ec(x, (o) can be rigorously defined in terms of an electromagnetic scattering matrix S for 
the cavity. This matrix relates incoming waves at wavevector k (frequency o) = ck) to 
outgoing waves in all of the asymptotic scattering channels and can be calculated from 
the wave equation. Note that while we speak of the frequency of the incoming wave in 
fact the S -matrix is a time-independent quantity depending on the wavevector k. This 
is the wavevector outside the cavity; in random lasers we will be interested in spatially 
varying dielectric functions so that in the "cavity" there is no single wavevector of the 
field. For any laser, including the random laser, the cavity can be defined as simply the 
surface of "last scattering" beyond which no backscattering occurs. The QB states are 
then the eigenvectors of the passive cavity S -matrix with eigenvalue equal to infinity; i. 
e. one has outgoing waves with no incoming waves. Because this boundary condition is 
incompatible with current conservation these eigenvectors have complex wavevector, 
k^i', these complex frequencies are the poles of the S -matrix and their imaginary parts 
must always be negative to satisfy causality conditions. There are normally a countably 
infinite set of such QB states. Due to their complex wavevector, asymptotically the 
QB states vary as r^('^^')/^exp(+|/m[^^]|r) and diverge at infinity, so they are not 
normalizable solutions of the time-independent wave equation. Therefore we see that 
QB states cannot represent the lasing modes of the cavity, even at threshold, as the 
lasing modes have real frequency and wavevector outside the cavity with conserved 
photon flux. 

When gain is added to the cavity the effect is to add another contribution to the 
dielectric function £g(x, ft)) which in general has a real and imaginary part. The imag- 
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FIGURE 14. Shift of the poles of the S-matrix in the complex plane onto 
the real axis to form threshold lasing modes when the imaginary part of 
the dielectric function e = Ec + £g varies for a simple Id edge-emitting cav- 
ity laser [33]. The cavity is a region of length L and uniform index n^- = 
1.5(fl), 1.05(i>)(ec = 2.25, 1.0025) terminated in vacuum at both ends. The cal- 
culations are based on the MB model discussed in section 5, with parameters 
k^L = 39 and Y± = 2. (a) = 1-5, the squares of different colors represent 
lm[e^,] = 0,-0.032,-0.064,-0.096,-0.128 (b) = 1.05,the squares of differ- 
ent colors represent Im[eg] = 0, —0.04, —0.08, —0.12, —0.16, note the increase in 
the frequency shift in the complex plane for the leakier cavity. The center of the 
gain curve is at kL = 39 which determines the visible line-pulling effect. 



inary pait of Eg has amplifying sign when the gain medium is inverted and depends 
on the pump strength; it compensates the outcoupling loss as well as any cavity loss 
from the cavity dielectric function e^- The specific form of this function for the MB 
model will be given in section 5 below. The threshold lasing modes (TLMs) are the so- 
lutions of the wave equation with Stotaii'^) = £c(x) + %(x) with only outgoing waves of 
real wavevector (we neglect henceforth for simplicity the frequency-dependence of 
ec(x)). The are the wavevectors of the TLMs with real lasing frequencies = ck^. 
These lasing wavevectors are clearly different from the complex ; moreover they are 
not equal to Re[k^] as often supposed. This can be seen by the following continuity ar- 
gument. Assume that £c(x) is purely real for simplicity, so that the S-matrix is unitary 
and all of its poles are complex and lie in the negative half plane. Turn on the pump, 
which we will call Do anticipating our later notation, so that the inversion rises steadily 
from zero, continuously increasing the amplifying part of Eg. The S-matrix is no longer 
unitary, and its poles move continuously "upward" towards the real axis until each of 
them crosses the axis at a particular pump value. Do (see Fig. 14); the place where each 
pole crosses is the real lasing frequency fc^ for that particular TLM. Note that the poles 
do not move vertically to reach the real axis but always have some shift of the lasing 
frequency from the passive cavity frequency, mainly due to line-pulling towards the 
gain center. As the Q-value of the cavity increases, the distance the poles need to move 
to reach the real axis decreases so that the frequency shift from can become 

very small and the conventional picture becomes more correct. In general the poles of 
the S-matrix are conserved quantities even in the presence of loss, so that the TLMs 
are in one to one correspondence with the QB states and thus are countably infinite, 
but for any cavity the pole which reaches the real axis first (i.e. at lowest pump Do) is 



the actual first lasing mode. At higher pump values the non-lineai^ effects of saturation 
and mode competition will affect the behavior; so only the lowest threshold TLM de- 
scribes an observable lasing mode for fixed pumping conditions, the first lasing mode 
at threshold. Which pole gets there first depends not only on the Q of the passive cav- 
ity resonance before gain is added, but also on the pai^ameters of £g(x) which include 
the atomic transition frequency, the gain linewidth and the pump conditions as will be 
discussed below. 

5. Self-consistent time-independent approachi to random lasing 

In section 4 we gave a general argument based on the scattering matrix with the addi- 
tion of gain to show that in general the QB states (passive cavity resonances) are never 
exactly the same as the threshold lasing modes (TLMs), even inside the cavity. How- 
ever the same argument indicated that inside a high-Q cavity the two sets of functions 
become very similar since the poles of the S -matrix are very close to the real axis and 
only a small amount of gain is required to move them to the real axis, which maps 
QB states onto TLMs. For localized states in the center of the sample the Q-values 
should be exponentially large and, as found numerically, QBs and TLMs should be 
indistinguishable (again, inside the cavity, outside the QB states have an unphysical 
growth). As already noted, the set of TLMs only defines threshold modes, as soon as 
the first TLM has turned on it will alter the gain medium for the other potential modes 
through spatial hole-burning and a non-linear approach needs to be considered. Very 
recently such an approach has been developed which has the major advantage of be- 
ing time-independent and partially analytic, providing both ease of computation and 
greater physical insight. The approach, due to Tiireci-Stone-Ge, is known as Ab Initio 
Self-consistent Laser Theory (AISC laser theory) [30,33,34]. It finds the stationary so- 
lutions of the MB semiclassical lasing equations in the multimode regime, for cavities 
of arbitrary complexity and openness, and to infinite order in the non-linear interac- 
tions. As such it is ideal for treating diffusive or quasi-ballistic random lasers which 
are extremely open and typically highly multimode even slightly above threshold. In 
this section we present the basic ideas with emphasis on threshold lasing modes, which 
is the focus of this review. The non-linear theory has been reviewed in some detail 
elsewhere [34], and we just present a brief introduction to it here. 

5. 1 . Maxwell-Bloch threshold lasing modes 

The MB semiclassical laser equations describe a gain medium of identical two-level 
"atoms" with energy level spacing hcOa = ficka and relaxation rate y\\ , being pumped by 
an external energy source. Do (which can vary in space), contained in a cavity which 
can be described by a linear dielectric function, £c(x). This leads to a population inver- 
sion of the atoms, D(x,t) which in the presence of an electric field creates a non-linear 
polarization of the atomic medium, P{x,t), which itself is coupled non-linearly to the 
inversion through the electric field, E{x,t). The electric field and the non-linear polar- 
ization are related linearly through Maxwell's wave equation, although above the first 
lasing threshold the polarization is implicitly a non-linear function of the electric field. 



The induced polarization also relaxes at a rate y±_ which is typically much greater than 
the rate at which the inversion relaxes, and this is a key assumption in our treatment 
of the non-linear regime, but will not be needed in the initial discussion of TLMs. 

The resulting system of non-lineai" coupled partial differential equations for the three 
fields E{x,t),P{x,t),D{x,t) are (c = 1): 

1 -V2£+-^P+ (10) 



£,.(x) e^x) 

p+ = -(i(o^ + Y^)p+ + ±-E+D (11) 

in 

D = (Do -D)- I (£+)*). (12) 

Here g is the dipole matrix element of the atoms and the units for the pump are cho- 
sen so that Do is equal to the time-independent inversion of the atomic system in the 
absence of an electric field. This pump can be non-uniform: Do = Do{x) based on the 
experimental pump conditions, but we will not discuss that case here. The electric field, 
polarization and inversion are real functions {E,P are vector functions in general, but 
we assume a geometry where they can be treated as scalars). In writing the equations 
above we have written these fields in the usual manner in terms of their positive and 
negative frequency components, E = £"+ + E^, P = P^ + P^, and then made the ro- 
tating wave approximation (RWA) in which the coupling of negative to positive com- 
ponents is neglected. There is no advantage in our treatment to making the standard 
slowly-varying envelope approximation and we do not make it. 

5.2. Self-consistent steady-state lasing equations 

The starting point of our formulation is to assume that there exists a steady state multi- 
periodic solution of equations (10)-(12) above, i. e., we try a solution of the form: 

E+{x,t)=f^^',ix)e-"^', P+{x,t)=f^P,{x)e-''^'. (13) 

Having taken c = 1 we do not distinguish between frequency and wavevector. The 
functions (x) are the unknown lasing modes and the real numbers are the un- 
known lasing frequencies; these functions and frequencies are not assumed to have 
any simple relationship to the QB states of the passive cavity and will be determined 
self-consistently. As the pump increases from zero the number of terms in the sum will 
vary, A^ = 0,l,2,... at a series of thresholds each new mode will appear. The general 
non-lineai- theory is based on a self-consistent equation which determines how many 
modes there are at a given pump, and solves for these modes and their frequencies. 
However in this section we will discuss threshold lasing modes (TLMs) and so we 
need only consider one term in the sum. Furthermore, at the first threshold the electric 
field is negligibly small and so the inversion is equal to the external pump profile, as- 
sumed uniform in space, D{x,t) = Dq. Assuming single-mode lasing the equation for 



the polarization becomes: 



-iDog^'¥a(x) 
fi[,y± - i[k^ - ka)) 

Having found (x) in terms of *Fp (x) , Dq , we substitute this result into the right hand 
side of Maxwell's equation along with (x) for the electric field on the left hand side. 
The result is: 

[V^ + Mx)iJ|>P„W = j^^^^ji^, (15) 

which can be written in the form: 

[V2 + (£,(x) + e,{x))kl]^>^{x) = 0, (16) 

where £g(x) is the dielectric function of the gain medium, which only varies in space 
if the external pump or the gain atoms are non-uniform. Defining convenient units of 
pump Dqc = fiy^/^Kk^g^ and replacing Do ^ Dq/Dqc, we find that 

e(x) = ^\ + ZM 1 (17) 

' krf^ + {k, - kaY + ri + - 

Equation (16) is to be solved with the boundary condition that at infinity one has 
only an outgoing wave at frequency k^, i. e., V,-*F^(x) = +/^^*F^(x) when r — ^ oo. In 
general this equation with this boundary condition cannot be solved for arbitrary choice 
of the lasing frequency k^ and for arbitrary values of the pump Dq; it is necessary to 

vary k^ and the pump strength Do to find the countably infinite set of values {k^,D\^^) 
at which a solution exists. This variation is equivalent to the pulling of the S-matrix 
poles onto the real axis discussed in section 4 above; d\^^ defines the threshold pump, 
for that pole, and k^ the point at which it crosses the real axis. As noted, while all 
of these solutions can be classified as TLMs, only the solution with the lowest value 
of Dq^^ will actually be a physical lasing state, as higher lasing modes are altered by 
non-linear modal interactions. 

Equation (16) shows that the TLMs are the solutions of the original Maxwell equa- 
tion with the addition of a complex, pump and frequency dependent dielectric function 
which is uniform in space (for the assumed uniform pumping). The imaginary and real 
parts of the gain dielectric function have the familiar symmetric and anti-symmetric 
two-level resonance form respectively. The dependence on the atomic frequency ka 
encodes the usual atomic line-pulling effect. In the limit of a very broad gain curve 
(y^ — )■ oo) the line -pulling effects can be neglected and we find the simple result 

Eg^-iDo/kl (18) 

i.e. a constant imaginary (amphfying) part of £^ proportional to the pump strength. 
Such linear gain models have been studied before, although typically with a constant 
imaginary part of the index of refraction instead of a constant imaginary part of the 



dielectric function. Our results show that in order to reproduce the TLMs of the MB 
equations one needs to take 



n{x) = edx) + eg{Do,k^ - ka,Y±) (19) 
so that the pump changes both the real and imaginary parts of the index of refraction. 
5.3. Solution for TLMs and CF states 
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FIGURE 15. Typical values of the threshold matrix elements ^^^^ in a two- 
dimensional random laser schematized in the inset Fig. 18 using sixteen CF states. 
The off-diagonal elements are one to two orders of magnitude smaller than the di- 
agonal ones. 



The differential equation (16) is self-consistent in the sense that the boundary condi- 
tions depend on the eigenvalue that one is solving for and so some form of non-linear 
search is required. The required search turns out to be much more convenient if one 
writes an equivalent integral form of the equation transforming it into a self-consistent 
eigenvalue problem. For this purpose we rewrite it in the form 



— p 

e,(x) 



(20) 



and then, treating the right hand side as a source, invert the equation with the appropri- 
ate Green function to obtain: 



iDor± ^fi f ^^Mx,x';k^)^>^{x') 



Y± - i{k^ - ka) kl Joj 



e,{x') 



(21) 



Here the integral is over the gain region which we will assume coincides with the cavity 
region The appropriate Green function satisfies 



[ec{x)-^V^ + k^]G{x,x'\k) = d'^{x-x'), 



(22) 




FIGURE 16. (a) False color plot of one TLM in a 2D random laser modeled as 
an aggregate of sub-wavelength particles of index of refraction n = 1.2 and ra- 
dius r ~ R/3Q against a background index n = 1 imbedded in a uniform disk of 
gain material of radius R(see inset. Panel (d)). The frequency of the lasing mode 
is kR = 59.9432, which is pulled from the real part of the dominating CF state 
ki„R = 59.8766 — 0.8593; (b) towards the transition frequency kaR ~ 60. The spa- 
tial profile of the TLM and CF state agree very well, whereas the corresponding 
QB state k„,R = 59.8602 - 0.8660/ (c) differs from that of the TLM and the CF 
state noticeably, as can be seen in Panel (d) where we plot the internal intensity 
along the 9 = 200° direction (white line in panel (a)). 



and is non-hermitian due to the outgoing wave boundary conditions: 
VrG{x,x'\k)\,-^co = V,-/G(x,x'|^)|,./^oo = ikG{x,x'\k), where V,. is the radial derivative. 
It has the spectral representation: 

C(x,xm=£54Jl^. (23) 

We refer to the functions (p,n{x,k) in (23) as the constant-flux (CF) states. They satisfy 

[£,{x)-'V^ + kl](pUx,k)=0 (24) 



with the corresponding non-hermitian boundary condition of purely outgoing spherical 
waves of fixed frequency k (eventually set equal to the lasing frequency) at infinity. 
Their dual (biorthogonal) partners (p„j{x',k) satisfy the complex conjugate differential 
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FIGURE 17. (a) False color plot of one TLM in a 2D random laser similar to 
that in the Fig. 16 but with particles of radius r = R/60 corresponding to weaker 
scattering (see inset. Panel (d)). The frequency of the lasing mode is kR = 29.9959, 
which is very close to the CF state k,„R = 30.0058 - 1 .3219? (b) but shifted from 
the coiTesponding QB state k„,R = 29.8813 - 1.3790/ (c). In Panel (d) we show 
the internal intensity of the three states in the 9 — n direction (white line in panel 
(a)); due to weaker scattering the QB state now differs substantially from the CF 
and TLM, which still agree quite well with each other 

equation with purely incoming wave boundary conditions. These dual sets satisfy the 
biorthogonahty relation: 

/ dx (pm (x, k) (p* (x,k) = 5,„n (25) 

with appropriate normalization. 

The CF states satisfy the standard wave equation Eq. (24), but with the non-hermitian 
boundary condition already mentioned; hence their eigenvalues kf^^ are complex, with 
(it can be shown) negative imaginary part, corresponding to amplification within the 
cavity. However outside the cavity, by construction, they have the real wavevector k^ 
and a conserved photon flux. They are a complete basis set for each lasing frequency k^ 
and hence they are a natural choice to represent the TLMs, as well as the lasing modes 



above threshold. Hence we make the expansion 



»P^(x)= £a^<(x) (26) 

m=l 

Substituting this expansion into Eq. (21), using biorthogonality, and truncating the ex- 
pansion at A'^ terms leads to the eigenvalue problem: 

Kt = DoAm [k^ ) / dx '— = ^^0 2^ Xipa>i , (27) 

where A,„{k) = iY^{k^ /kl)/\{y^ - i{k - k,)){e - ki{k))]. 

One sees that the TLMs in the CF basis are determined by the condition that an 
eigenvalue of the matrix Do3^^^\k^) is equal to unity. Since the matrix 3^^^\k^) is 
independent of Do it is natural to focus on this object, which we call the threshold 
matrix. It is a complex matrix with no special symmetries, implying that its eigenvalues 
are all complex for a general value of k^. If the real control parameter Do (the 
pump) is set equal to 1/|A^| then the matrix Dq,'^^^^'' {kp) will have an eigenvalue of 
modulus unity, but not a real eigenvalue equal to unity as required and no solution 
for the TLMs exists for this choice of k^ . It is the phase condition, that (^^ ) must 
be real that determines the allowed lasing frequencies. In practice one orders the 
in decreasing modulus based on an initial approximation to the lasing frequency, k^, 
and then tunes k^ slowly until each eigenvalue flows through the real axis (which is 
guaranteed by the dominant k-dependence contained in the factor A„,(^)). Normally 
the eigenvalues do not switch order during this flow and the largest eigenvalue 

will determine the lowest threshold TLM, with threshold D^^^ = l/X^{k^), where k^ 
is the frequency which makes the largest eigenvalue J^^^\k^) real. The eigenvector 
corresponding to A^ gives the coefficients for the CF expansion of the TLM of the 
first mode ^'^{x). TLMs with higher thresholds can be found by imposing the reality 
condition on smaller eigenvalues of ^^^\k^). This approach has been described in 
detail elsewhere [32,34], and provides a much more efficient method for finding TLMs 
than solving the self-consistent differential equation, (16). 

We immediately see from Eqs. (25) and (27) that for an arbitrarily shaped cavity of 
uniform dielectric constant Sc the matrix 2^^^\k^) is diagonal due to biorthogonality 
of the CF states. Thus each TLM is a single CF state, coiTcsponding to one of the 
k^i which satisfies the reality condition. In this case the expansion of (x) consists 
of just one term and the threshold lasing equation is equivalent to the Eq. 24 with 
appropriate relabeling. When varies in space, as for random lasers, the threshold 
matrix is not diagonal and there can in principle be many CF states contributing to one 
TLM. However since (p„, (x) , (x) are uncorrelated fluctuating functions of space, it 
turns out that the threshold matrix in RLs is approximately diagonal and the threshold 
modes are dominated by one, pseudo-random CF state determined by solving Eq. 24 
for the appropriate random dielectric function ec{x). This is shown in Fig. 15 below. 
In summary, the theory leading to the threshold equation (27) gives an efficient time- 
independent method for finding the TLMs of random lasers in any disorder regime. In 



general these TLMs are very close to a single CF state determined by the Eq. 24 at the 
lasing frequency ^^.With this new method TLMs of random lasers can be found for 
complex two and even three-dimensional geometries. In Figs. 16 and 17 we compare 
TLMs, CF states and QB states for the two-dimensional random laser model used in 
ref. [32], illustrating the agreement of TLMs with CF states even for weak scattering, 
while a significant deviation from the closest QB state is found. 

This AISC laser theory is well-suited to describe not just TLMs but to find the true 
multimode lasing spectrum of RLs above threshold. This will not be treated in detail 
here, but in the next section we briefly explain the basic approach in the non-linear 
theory and show one representative result. 

5.4. Non-linear AISC laser theory 

The key to generalizing this theory to the multi-mode non-linear regime is to return 
to the fundamental MB equations and go beyond the assumption that the inversion 
D{x,t) is equal to the constant threshold pump Dq. Once lasing modes have turned 
on their spatially varying electric fields cause varying degrees of stimulated emission 
from the gain atoms and hence tends to reduce the inversion D from the pump value 
Dq in a manner which varies in space and in principle in time. However it has been 
shown that if Y\\, then the time-dependence of the inversion is weak and although 
D is varies in space, it is a good approximation to take D{x,t) = D{x). This stationary 
inversion approximation (SIA) has been used in laser theory for many years, going back 
to Haken [60], but has not been incorporated into an ab initio method such as AISC 
laser theory. We will not review the details of the derivation of the non-linear multimode 
theory of Tiireci-Stone-Ge, which have been given elsewhere [30, 34]. Instead we just 
state that the net effect of the non-linear interactions within the SIA is just to replace 
the uniform inversion as follows: 

in all of the equations of the theory of the TLMs. Here v labels all above thi^eshold 
modes and r{ky) is a Lorentzian centered at the lasing frequency of mode v with width 
Y±. If we make this substitution into the Eq. (21) we arrive at the fundamental integral 
equation of AISC laser theory: 

^ fx) = -JR^I^i f , . Gix,x'-k,)^',ix') 

^ y^-i{h-ka)kl.loj £,(x')(l+Ivrv|»I'v(x')P)' ^ ' 

Note that this equation shows that each lasing modes interacts with itself (saturation) 
and all other lasing modes (mode competition) via the "hole-burning" denominator of 
Eq. (28). This set of coupled non-linear equations is still conveniently solved in the 
basis of CF states for each modal frequency as for the TLMs; the details have been 
given elsewhere [32,34]. 

The first results of the AISC laser theory for the modal properties of multimode 
random lasers in weak-scattering two-dimensional media were given in [32]. We will 
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FIGURE 18. (a) CF (dots) and QB (crosses) frequencies in a 2D random laser 
modeled as an aggregate of sub-wavelength particles of index of refraction n = 1 .2 
against a background index n = I imbedded in a uniform disk of gain material 
(see inset). The two sets of complex frequencies are statistically similar but dif- 
fer substantially. The solid curve shows the gain curve r(^) with 7x = 1. (b) 
Lasing frequencies of the same random system well above threshold (coloured 
lines). Coloured circles denote the CF state dominating the correspondingly col- 
ored modes at threshold. 

not present a full picture of these results here, but just show some properties of the RL 
lasing frequencies in Fig. 1 8 below. The model is explained in the figure caption (see 
inset). The complex CF and QB frequencies are shown to be distinct and the lasing 
frequencies are subject to very strong line-pulling effects. 

The new tool of AISC laser theory allows one to study random lasers with full non- 
linear interactions in 2D and even in 3D. The elimination of time-dependence in this 
theory makes larger and more complex cavities computationally tractable. The theory 
also provides a new language based on CF states to describe the lasing modes. Now 
detailed statistical studies as well as comparisons to statistical models based on random 
matrix theory, disordered media theory and wave chaos theory are needed. Such studies 
are in progress. 

6. Conclusion 

A decade of theoretical study of random lasers has clarified the nature of the lasing 
modes in disordered systems with multiple scattering and gain. Most importantly it 
has been established that high-Q passive cavity modes such as those created by Ander- 
son Localization or by rare fluctuations of various kinds are not necessary in order to 
have self-organized laser oscillation at a frequency distinct from the atomic transition 
frequency (gain center). In addition this study has emphasized a point of general impor- 
tance in laser theory, that threshold lasing modes are not identical to the quasi-bound 
states (resonances) of the passive cavity. This point is demonstrated by a number of nu- 
merical calculations presented above and also can be understood from the realization 
that the QB states are eigenvectors of the unitary S-matrix of the cavity without gain, 
but at complex frequency, whereas the threshold lasing modes are eigenvectors of the 
non-unitary S-matrix of the cavity with gain and with real frequency. The difference 
between these eigenvectors (within the cavity), which is large in the weak scattering 



limit, becomes small in the diffusive regime as the Q of the cavity increases and is 
negligible, e.g. for Anderson localized modes and for high-Q modes of conventional 
cavities. The new basis set of constant flux states provides a better approximation for 
finding the threshold lasing modes of random lasers and coincides with the exact lasing 
modes of uniform index cavities. Further statistical and analytical study is necessary to 
characterize the properties of random lasers in the different regimes, weak scattering, 
diffusive and localized, and to understand the effects of non-linear interactions. 
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A. The multipole method 

This appendix details the principle of the multipole method as used in this paper and its 
implementation. Although we describe here the method for two-dimensional systems, 
it can be also applied to three-dimensional structures. 

We consider a random collection of Nc non-overlapping cylinders with arbitrary 
complex dielectric constant £/ = £/' + /£/" = nj and arbitrary radii ai located in a uni- 
form medium with complex dielectric constant Sb = £b + i£b" = n\ (Fig. 19), where 
«/ = «J + in^ and = + are the refractive indices of the cyhnders and the back- 
ground. The complex dielectric permittivities of the cylinders and the background can 
be arbitrary and may be frequency dependent. 

In two dimensions, the solution of the electromagnetic field problem decouples into 
two fundamental polarizations, in each of which the field may be characterized by a 
single field component: V{y) = (for TM polarization) and V{r) = H, (for TE polar- 
ization). In the co-ordinate system that is used, the z axis is aligned with the cylinder 
axes. 

The field component V satisfies the Helmholtz equation 

V2y(r) + A:V(r)V(r) =0. (30) 

For TM polarization, both V{r) and its normal derivative v • W(r) are continuous 
across all boundaries, while for TE polarization the corresponding boundary conditions 
are the continuity of V(r) and its weighted normal derivative v • W/n^(r). Here, n{r) 
denotes the refractive index of the relative medium and v is an unit outward normal 
vector. 

In the vicinity of the /"^ cylinder, we may represent the exterior field in the back- 
ground medium (refractive index ni,) in local coordinates ri = {ri,di) = r — ci where C/ 



o 



FIGURE 19. The geometry and local co-ordinate systems. 



represents the center of the cylinder and write 



V(r)= £ \A[Mknhr,)+BlHl^\knhri)\e" 



JmBi 



(31) 



This local expansion is valid only in an annulus extending from the surface of the 
cylinder / to the surface of the nearest adjacent cylinder. 

The global field expansion (also refeiTed to as a Wijngaard expansion), which is valid 
everywhere in the background matrix, comprises only outgoing cylindrical harmonic 
terms: 



Correspondingly, the field inside any cylinder / is written in an interior expansion: 



Then, applying Graf's addition theorem [50] to the terms on the right hand side of 
(32) (see Fig. 19), we may express the global field expansion in terms of the local 
coordinate system for the cylinder. Equating this with the local expansion (31), we 
deduce the field identity (also known as the Rayleigh identity): 




(32) 



q=l in= 




(33) 



A 



m 



(34) 



where 




(35) 



Here, [ciq^diq) are the polar coordinates of the vector c;^ = — c/, the position of 
cylinder q relative to cylinder I. 

This is the first connection between the standing wave ({A,'„}) and outgoing ({5^,,}) 
multipole coefficients, one which follows solely from the system geometry. Eq. (34) 
indicates that the local field in the vicinity of cylinder / is due to sources on all other 
cylinders {q ^ I), the contributions of which to the multipole term of order m — p at 
cylinder / are given by H^p. 

The second relation between the {A'^} and multipole coefficients is obtained 
from the field continuity equations (i.e., the boundary conditions) at the interface of 
cylinder / and the local exterior (31) and interior field (33) expansions. From these, we 
obtain: 

B'm = R'm^L^ (36) 

where the interface reflection and transmission coefficients, for both and polar- 
ization are given by 

; _ ^niJ^„{kniai)J,niknhai) - ni,Jni{kniai)J',„{knbai) 

f (I) (1)' ' V-'"/ 

c,njJ'^{kniai)H,n ' {kribai) - nhJm{kniai)H,n {hihai) 

T' = 2i/{nkaL) 

C,niJ'„{nikai)Hm ' {kuhUi) - ni,Jm{kniai)Hm {kubOi) 



in which = 1 for TM polarization and ^ = ?i^(r)/?i^(r) for TE polarization. 

To derive a simple closed form expression for the solution of the problem, we use 
partitioned matrix notation, introducing vectors A' = and B' = [B^^] and expressing 
(34) in the form 

A'=£h''?B^ (40) 

<? 

where A' and B' denote vectors of multipole coefficients for cylinder /. The matrix H 
is block partitioned according to H'^ = [H^p] for I ^ q (35), and H" = [0], each block 
of which is a matrix of Toeplitz form. CoiTcspondingly, the matrix forms of Equations 
(36) and (37) are 

B = RA, (41) 
C = TA, (42) 

where R = diagR' is a block diagonal matrix of diagonal matrices R' = diag/?^, and 
with corresponding definitions applying for the transmission matrices. 

Then, with the introduction of the partitioned vectors A = [A'], B = [B'] and the 
partitioned matrix H = [H''^] , we form the system of equations 

(I-RH)B = 0. (43) 



The problem has now been reduced to the solution of a generalized eigenvalue prob- 
lem for the matrix equation (43). The nontrivial solutions of the secular equation (43) 



determines modes of the random system. Finding the nontrivial solutions of the hnear 
system of equations (43) requires that the determinant of the system matrix vanishes: 
(44) 

D = 0, where D = det(S"^) with (44) 
S-i(A) = (I-RH). (45) 

Equivalently, this problem may be recast as a search for the poles of the scattering ma- 
trix S(A) (i.e., solutions of detS^^ (A) = 0). Once the pole is located, the con^esponding 
null vectors B of (43) are the multipole coefficients of the scattered field which are used 
to calculate the QB state profiles exterior to the scatterers using (32). The field inside 
a cylinder is calculated according to the interior expansions (42) and(33). The TLM 
poles must be searched in the (A,£c) domain, given the pump changes not only the 
imaginary part of the refractive index but the real part as well (4). 

The formal system (44) is of infinite dimension and so must be truncated in order to 
generate a computational solution, the accuracy of which is governed by the number of 
retained multipole coefficients N,n = INn^^a + 1 > where A'max is the truncation order of 
the multipole series, i.e., only the terms corresponding to the cylindrical harmonics of 
order n = —Nmax, ■■■,Nmax are retained. 

B. Finite Element IVIetliod 

We have also used the Finite Element Method (FEM) [57], implemented in a commer- 
cial software (ComsoF^)^ to solve the wave equation (30) and calculate the complex 
eigenvalues and eigenfunctions of the passive modes of the systems that were calcu- 
lated by the multipole method. The method suitably applies for modeling passive or 
active modes in a cavity, which is suiTounded by perfectly matched layers [61] to sim- 
ulate open boundaries. It is possible to obtain all the leaky modes, even the resonances 
characterized by a very small quality factor (as small as 5), in a reasonable computation 
time with a commercial PC, provided the size of the geometry is smaller than hundred 
times the wavelength. This is in contrast with the other methods described in this paper, 
which require much heavier computation. 

One of the most important step of the Finite Element Method is the creation of 
the mesh which describes the system. Figure 20 shows a close up on a typical mesh 
calculated for the 2D random system of Fig. 2. The maximum size of elements must be 
smaller than 7 times the wavelength [62]. 




FIGURE 20. Close up on a typical mesh created by Comsor"^to describe the 2D 
random system of Fig. 2. 



References 

1. D. S. Wiersma, "The smallest random laser," Nature 406, 132 (2000). 

2. H. Cao, "Random Lasers with Coherent Feedback," in Optical Properties ofNanostructured Random 
Media V. M. Shalaev editor. Topics Appl. Phys. 82, 303330, (Springer- Verlag Berlin Heidelberg 

2002) . 

3. H. Cao, "Lasing in random media," Waves in Random Media 13, Rl (2003). 

4. H. Cao, J. Y. Xu, Y. Ling, A. L. Burin, E. W. Seeling, X. Liu, and R. R H. Chang, "Random Lasers 
With Coherent Feedback," IEEE J.Selected Topics in Quantum Elec, 9, No. 1, (January/February 

2003) . 

5. D. S. Wiersma, "The physics and applications of random lasers," Nature Physics 4, 359-367 (2008). 

6. P. W. Anderson, "Absence of Diffusion in Certain Random Lattices," Phys. Rev. 109, 1492 (1958). 

7. C. Gouedart, D. Husson, C. Sauteret, F. Auzel, and A. Migus, "Generation of spatially incoherent 
short pulses in laser pumped neodymium stoichiometric crystal powders," J. Opt. Soc. Am. B 10, 
2358 (1993). 

8. N. M. Lawandy, R. M. Balachandra, A. S. L. Gomes, and E. Sauvain, "Laser Action in strongly 
scattering media," Nature 368, 436 (March 1994). 

9. W. L. Sha, C. H. Liu, and R. R. Alfano, "Spectral and temporal measurements of laser action of 
Rhodamine 640 dye in strongly scattering medi". Opt. Lett. 19, 1922 (1994). 

10. H. Cao, Y.G. Zhao, H. C. Ong, S. T. Ho, J. Y. Dai, J. Y. Wu, and R. R H. Chang, "Ultraviolet lasing in 
resonators formed by scattering in semiconductor polycristalline films," Appl. Phys. Lett. 73, 3656 
(1998). 

11. H. Cao, Y. G. Zhao, S. T. Ho, E. W. Seelig, Q. H. Wang, and R. R H. Chang, "Random Laser Action 
in Semiconductor Powder," Phys. Rev. Lett. 82, 2278 (1999); H. Cao, Y. G. Zhao, H. C. Ong, and R. 
P. H. Chang, "Far-field characteristics of random lasers", Phys. Rev. B 59, 15107 (1999). 

12. S. V. Frolov, Z. V. Vardeny, K. Yoshino, A. Zakhidov, and R. H. Baughman, "Stimulated emission 
in high-gain organic media ," Phys. Rev B 59, 5284 (1999). 

13. H. Cao, J. Y. Xu, D.Z. Zhang, S.-H. Chang, S. T. Ho, E. W. Seehg, X. Liu and R. R H. Chang, 
"Spatial confinement of laser light in active random media," Phys. Rev. Lett. 84, 5584 (2000). 

14. H. Cao, Y. Ling, J. Y. Xu, C. Q. Cao and C. Q. Cao, "Photon Statistics of Random Lasers with 
resonant feedback", Phys. Rev. Lett. 86, 4524 (2001). 

15. A. Lagendijk, B. van Tiggelen, and D. S. Wiersma, "Fifty years of Anderson localization," Phys. 
Today 62, (2009) 24-29. 

16. V. Milner, and A. Z. Genack, Photon Localization Laser: "Low-Threshold Lasing in a Random 
Amplifying Layered Medium via Wave Localization," Phys. Rev. Lett. 94, 073901 (2005). 

17. K. L. van der Molen, R. W. Tjerkstra, A. P. Mosk, and A. Lagendijk, "Spatial Extent of Random 
Laser Modes," Phys. Rev. Lett. 98, 143901 (2007). 

18. V. M. Apalkov, M. E. Raikh and B. Shapiro, "Random Resonators and Prelocalized Modes in Dis- 
ordered Dielectric Films," Phys. Rev. Lett. 89, 016802 (2002). 

19. V. M. Apalkov and M. E. Raikh, and B. Shapiro, "Almost localized photon modes in continuous and 
discrete models of disordered media", J. Opt. Soc. Am. B Vol. 21, (2004). 

20. S. Mujumdar, M. Ricci, R. Torre, D. S. Wiersma, "Amplified Extended Modes in Random Lasers," 
Phys. Rev. Lett. 93, 053903 (2004). 

21. S. Mujumdar, V. Trek, R. Torre, and D. S. Wiersma, "Chaotic behavior of a random laser with static 
disorder," Phys. Rev. A76, 033807 (2007). 

22. A. A. Chabanov, Z.Q. Zhang, and A. Z. Genack, "Breakdown of Diffusion in Dynamics of Extend- 
edWaves in Mesoscopic Media," Phys. Rev. Lett. 90, 203903 (2003). 

23. J. Fallert, R. J. B. Dietz, J. Sartor, D. Schneider, C. Klingshim, and H. Kalt, "Co-existence of strongly 
and weakly localized random laser modes". Nature Photon. 3, 279282 (2009). 

24. X. Jiang and CM. Soukoulis, "Time dependent Theory for Random Lasers," Phys. Rev. Lett. 85, 70 
(2000). 



25. C. Vanneste, and P. Sebbah, "Selective Excitation of Localized Modes in Active Random Media," 
Phys. Rev. Lett. 87, 183903 (2001). 

26. P. Sebbah and C. Vanneste, "Random laser in the localized regime," Phys. Rev. B66, 144202 (2002). 

27. C. Vanneste, P. Sebbah, and H. Cao, "Lasing with Resonant Feedback inWeakly Scattering Random 
Systems," Phys. Rev. Lett. 98, 143902 (2007). 

28. C. Vanneste, and P. Sebbah, "Complexity of two-dimensional quasimodes at the transition from 
weak scattering to Anderson localization", Phys. Rev. A79, 041802(R) (2009). 

29. X. Wu, W. Fang, A. Yamilov, A. A. Chabanov, A. A. Asatryan, L. C. Botten, and H. Cao, "Random 
lasing in weakly scattering systems," Phys. Rev. A74, 053812 (2006). 

30. H. E. Tureci, A. D. Stone, B. Collier, "Self-consistent multimode lasing theory for complex or ran- 
dom lasing media," Phys. Rev. A 74, 043822 (2006). 

31. H. E. Tureci, A. D. Stone, L. Ge, "Theory of the spatial structure of nonlinear lasing modes," Phys. 
Rev. A 76, 013813 (2007). 

32. H. E. Tureci, L. Ge, S. Rotter, A. D. Stone, "Strong Interactions in Multimode Random Lasers," 
Science 320, 643 (2008). 

33. L. Ge, R. Tandy, A. D. Stone, and H. E. Tureci, "Quantitative verification of ab initio self-consistent 
laser theory," Optics Express, 16, 16895 (2008). 

34. H. E. Tureci, A. D. Stone, L. Ge, S. Rotter and R. J. Tandy, Non-Linearity, "Afo initio self-consistent 
laser theory and random lasers," 22, CI -CI 8 (2009). 

35. P. Sebbah, D. Somette and C. Vanneste, "Anomalous diffusion in two-dimensional Anderson- 
localization dynamics," Phys.Rev. B48, 12506 (1993). 

36. F. A. Pinheiro, M. Rusek, A. Orlowski, and B. A. van Tiggelen, "Probing Anderson localization of 
light via decay rate statistic," Phys.Rev. E 69, 026605 (2004). 

37. The idea of considering localized modes for random lasing can already be found in P. Sebbah, D. 
Somette, and C. Vanneste, OSA Proc. on Advances in Optical Imaging and Photon Migration, ed. 
R.R. Alfano, Vol. 21, 68 (1994) and P. Sebbah, "A new approach for the study of wave propagation 
and localization," PhD. Thesis, Universite de Nice - Sophia Antipolis, France (1993). 

38. A. E. Siegman, Lasers (University Science Books, Mill Valley, 1986). 

39. A. Taflove, Computational Electrodynamics: The Finite-Difference Time-Domain Method (Artech 
House, Norwood, 1995). 

40. J. P. Berenger, "A perfectly matched layer for the absorption of electromagnetic waves," J. Comput. 
Phys. 114, 185-200(1995). 

41. X. Jiang and C. M. Soukoulis, "Localized random lasing modes and a path for observing localiza- 
tion," Phys. Rev. E 65, 025601(R) (2002). 

42. A. G. Fox and T. Li, "Resonant modes in an optical maser," Proc. IRE 48, 1904-1905 (1960); "Res- 
onant modes in a maser interferometer," Bell Syst. Tech. J. 40 453-488 (1961). 

43. S. M. Dutra and G. Nienhuis, "Quantized mode of a leaky cavity," Phys. Rev. A62, 063805 (2000). 

44. J. C. J. Paasschens, T. Sh. Misirpashaev, and C. W. J. Beenakker, "Localization of light : Dual 
symmetry between absorption and amplification," Phys. Rev. B 54, 1 1887 (1996). 

45. A. A. Asatryan, N. A. Nicorovici, L. C. Botten, M. C. de Sterke, P. A. Robinson, and R. C. McPhe- 
dran, "Electromagnetic localization in dispersive stratified media with random loss and gain," Phys- 
ical Review B, 57, 13535-13549 (1998). 

46. X. Jiang, Q. Li, and CM. Soukoulis, "Symmetry between absorption and amplification in disordered 
media," Phys. Rev. B 59, R9007 (1999-11). 

47. X. Wu, J. Andreasen, H. Cao, and A. Yamilov, "Effect of local pumping on random laser modes in 
one dimension," J. Opt. Soc. Am. B 24, A26-A33 (2007). 

48. D. Felbacq, G. Tayeb, D. Maystre, "Scattering by a random set of parallel cylinders," J. Opt. Soc. 
Am. A, 11, pp.2526-2538 (1994). 

49. E. Centeno, and D. Felbacq, "Characterization of defect modes in finite bidimensional photonic 
crystals," J. Opt. Soc. Am A 16, pp.2705-2712 (1999). 

50. A. A. Asatryan, K. Busch, R. C. McPhedran, L. C. Botten, C. M. de Sterke, and N. A. Nicorovici, 
"Two-dimensional Green tensor and local density of states in finite-sized two-dimensional photonic 
crystals," Waves Random Media 13, 9-25 (2003). 

51. A. A. Asatryan, K. Busch, R. C. McPhedran, L. C. Botten, C. Martijn de Sterke, and N. A. 



Nicorovici, "Two-dimensional Greens function and local density of states in photonic crystals con- 
sisting of a finite number of cylinders of infinite length," Phys. Rev. E, 63, 046612 (2001). 

52. A. A. Asatryan, L. C. Botten, N. A. Nicorovici, R. C. McPhedran, C. M. de Sterke, "Frequency shift 
of sources embedded in finite two-dimensional photonic clusters," Waves in Random and Complex 
Media, 16, 151-165 (2006). 

53. K. M. Lo, R. C. McPhedran, I. M. Bassett and G. W. Milton, "An electromagnetic theory of optical 
wave-guides with multiple embedded cylinders," J. of Lightwave Tech., 12, 396-410 (1994). 

54. T. P. White, B. Kuhlmey, R. C. McPhedran, D. Maystre, G. Renversez, C. Martijn de Sterke, and 
L.C. Botten, "Multipole method for microstructured optical fibers I: Formulation" J. Opt. Soc. Am. 
B 10, 2322-2330 (2002). 

55. B. Kuhlmey, T. P. White, G. Renversez, D. Maystre, L.C. Botten, C. Martijn de Sterke, and R. C. 
McPhedran, "Multipole method for microstructured optical fibers IL Implementation and results," J. 
Opt. Soc. Am. B, 10, 2331-2340, (2002). 

56. D.P. Fussell, R.C. McPhedran, C. Martijn de Sterke, and A.A. Asatryan, "Three- dimensional local 
density of states in a finite two-dimensional photonic crystal composed of cylinders," Phys. Rev. E 
67, 045601(R) (2003). 

57. J. Jin, The finite element method in electromagnetics (Wiley, New York, 1993). 

58. EPAPS link to four movies representing the space-time variation of the QB states and TLM over a 
few cycles in two different regimes, localized {n' = 2.0) and diffusive (n' = 1.25). 

59. M. Sargent, M.O. Scully, and W.E. Lamb, Laser Physics (Addison- Wesley Publishing Company, 
Reading, 1982). 

60. H. Haken, Light: Laser Dynamics (North-Holland, Amsterdam, 1985), vol. 2. 

61. U. Pekel, R. Mittra, "A finite-element-method frequency-domain application of the perfectly 
matched layer (PML) concept," Microwave and optical technology letters Vol 9, Issue 3, 2007. 

62. "Modelling and characterization of air/silica microstructured fibres for optical telecommunication 
applications," PhD, Peyrilloux Ambre, 2003, Limoges (France) 



0.04 
0.035 

0.03 
0.025 

0.02 
0.015 

0.01 
0.005 




'MODE_D4558g.Y' using 2:3 
'MODE4567_0_0078.Y' using 2:3 



I 



